Visual Analytics exercise on Armed conflicts - Initial Draft

On this page

• 4.1 The Task
• 4.2 Getting Started
  • 4.2.1 Loading R packages
  • 4.2.2 Importing and loading the ACLED data
  • 4.2.3 Downloading and loading the shape files for country
• 4.3 Data Preparation and Cleaning
  • 4.3.1 Loading Admin 1 & 2 (administrative region/area) shape files
  • 4.3.2 Data Wrangling
• 4.4 Choropleth Maps
  • 4.4.1 Choropleth map of Incidents & Fatalities by Admin2 level (by District)
• 4.5 Proportional Symbol Maps
  • 4.5.1 Visualising the location of conflict events
• 4.6 UI design - Part 1
• 4.7 Data Preparation for Spatial Analysis
  • 4.7.1 Filtering the Event and Year (Event type = Battles, in 2023)
  • 4.7.2 Computing Contiguity Spatial Weights
  • 4.7.3 Row-standardised weights matrix
• 4.8 Cluster and Outlier Analysis
  • 4.8.1 Computing Local Moran’s I
  • 4.8.2 Mapping the Local Moran’s I
  • 4.8.3 Mapping Local Moran’s I values
  • 4.8.4 Mapping Local Moran’s I p-values
  • 4.8.6 Data table for the Moran’s I values
• 4.9 UI design - Part 2
• 4.10 Creating the LISA cluster map
  • 4.10.1 Plotting the Moran Scatter plot
  • 4.10.2 Plotting Moran scatterplot with standardised variable
  • 4.10.3 Preparing LISA map classes
  • 4.10.4 Plotting LISA Map
• 4.11 UI design - Part 3
• 4.12 Hot & Cold Spot Area Analysis
  • 4.12.1 Deriving distance-based weight matrix
    • 4.12.1.1 Deriving the centroid
    • 4.12.1.2 Determining the cut-off distance
  • 4.12.2 Computing Adaptive distance weight matrix
  • 4.12.3 Computing GI statistics - Fixed distance
    • 4.12.3.1 Mapping Gi values with Fixed distance weights
  • 4.12.4 Computing GI statistics - Adaptive distance
    • 4.12.4.1 Mapping Gi values with Adaptive distance weights
• 4.13 UI design - Part 4
• 4.14 References

4.1 The Task

In this take-home exercise, we are required to select one of the modules of our proposed Shiny application (Group Project) and complete the following tasks:

• To evaluate and determine the necessary R packages needed for our Shiny application are supported in R CRAN,

• To prepare and test that the specific R codes can run and returns the correct output as expected,

• To determine the parameters and outputs that will be exposed on the Shiny applications,

• To select the appropriate Shiny UI components for exposing the parameters determined above, and

• To include a section called UI design for the different components of the UIs for the proposed design.

4.2 Getting Started

Our project will be using open-source data from the Armed Conflict Location & Event Data Project (ACLED).

Specifically, our project will be focusing on the visualisation of Armed conflicts in Myanmar, and I will be preparing the modules on Cluster & Outlier Analysis (LISA) and Hot/Cold zone analysis.

4.2.1 Loading R packages

The below R packages will be used in this exercise and for the Shiny application

pacman::p_load(sf, tidyverse, tmap, dplyr,
               spatstat, spdep,
               lubridate, leaflet,
               plotly, DT, viridis,
               ggplot2, sfdep)

4.2.2 Importing and loading the ACLED data

Country specific data from the Armed Conflict Location & Event Data Project (ACLED) can be downloaded at https://acleddata.com/data-export-tool/

Loading the ACLED data set for Myanmar.

ACLED_MMR <- read_csv("data/MMR.csv")

4.2.3 Downloading and loading the shape files for country

Shape files were downloaded from the Myanmmar Information Management Unit (MIMU) website.

I chose this source over GADM and GeoBoundaries due to its updated administrative region information and map levels.

Note - Data Quality Issues

ACLED captures event data from national, sub-national and other media sources, and populates event locations based on the last known information.

However, some names of administrative areas were found to have changed; either disaggregated into new administrative areas or previously active but now defunct. Some administrative areas were also aggregated into higher administrative areas.

As part of our data cleaning and preparation process, I had to identify discrepancies in both admin1 & 2 (administrative levels) and re-name some administrative areas to sync with the downloaded shape files from MIMU.

4.3 Data Preparation and Cleaning

I will first load in the shape files at the admin 1 (region/state) and admin 2 (districts) levels. Most of our plots will be utilizing admin 2 levels.

4.3.1 Loading Admin 1 & 2 (administrative region/area) shape files

mmr_shp_mimu_1 <-  st_read(dsn = "data/geospatial3",  
                  layer = "mmr_polbnda2_adm1_250k_mimu_1")

Reading layer `mmr_polbnda2_adm1_250k_mimu_1' from data source 
  `C:\imranmi\imran's data sc\R-ex\R-Ex5\data\geospatial3' using driver `ESRI Shapefile'
Simple feature collection with 18 features and 6 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS:  WGS 84

mmr_shp_mimu_2 <-  st_read(dsn = "data/geospatial3",  
                  layer = "mmr_polbnda_adm2_250k_mimu")

Reading layer `mmr_polbnda_adm2_250k_mimu' from data source 
  `C:\imranmi\imran's data sc\R-ex\R-Ex5\data\geospatial3' using driver `ESRI Shapefile'
Simple feature collection with 80 features and 7 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS:  WGS 84

class(mmr_shp_mimu_2)

[1] "sf"         "data.frame"

The Shape file for admin2 level map, is an SF object, with geometry type: Multipolygon.

st_geometry(mmr_shp_mimu_2)

Geometry set for 80 features 
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS:  WGS 84
First 5 geometries:

st_crs(mmr_shp_mimu_2)

Coordinate Reference System:
  User input: WGS 84 
  wkt:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]],
        ID["EPSG",6326]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433],
        ID["EPSG",8901]],
    CS[ellipsoidal,2],
        AXIS["geodetic longitude",east,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic latitude",north,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]]]

Next, I will check the unique district names in this shape file (admin2)

unique_regions_mimu2 <- unique(mmr_shp_mimu_2$DT)

unique_regions_mimu2

 [1] "Hinthada"                        "Labutta"                        
 [3] "Maubin"                          "Myaungmya"                      
 [5] "Pathein"                         "Pyapon"                         
 [7] "Bago"                            "Taungoo"                        
 [9] "Pyay"                            "Thayarwady"                     
[11] "Falam"                           "Hakha"                          
[13] "Matupi"                          "Mindat"                         
[15] "Bhamo"                           "Mohnyin"                        
[17] "Myitkyina"                       "Puta-O"                         
[19] "Bawlake"                         "Loikaw"                         
[21] "Hpa-An"                          "Hpapun"                         
[23] "Kawkareik"                       "Myawaddy"                       
[25] "Gangaw"                          "Magway"                         
[27] "Minbu"                           "Pakokku"                        
[29] "Thayet"                          "Kyaukse"                        
[31] "Maungdaw"                        "Mrauk-U"                        
[33] "Sittwe"                          "Thandwe"                        
[35] "Hkamti"                          "Kale"                           
[37] "Kanbalu"                         "Katha"                          
[39] "Kawlin"                          "Mawlaik"                        
[41] "Monywa"                          "Naga Self-Administered Zone"    
[43] "Sagaing"                         "Shwebo"                         
[45] "Tamu"                            "Yinmarbin"                      
[47] "Kengtung"                        "Monghsat"                       
[49] "Tachileik"                       "Hopang"                         
[51] "Kokang Self-Administered Zone"   "Kyaukme"                        
[53] "Lashio"                          "Matman"                         
[55] "Mongmit"                         "Muse"                           
[57] "Pa Laung Self-Administered Zone" "Danu Self-Administered Zone"    
[59] "Langkho"                         "Loilen"                         
[61] "Pa-O Self-Administered Zone"     "Taunggyi"                       
[63] "Dawei"                           "Kawthoung"                      
[65] "Mandalay"                        "Meiktila"                       
[67] "Myingyan"                        "Nyaung-U"                       
[69] "Pyinoolwin"                      "Yamethin"                       
[71] "Mawlamyine"                      "Thaton"                         
[73] "Det Khi Na"                      "Oke Ta Ra"                      
[75] "Kyaukpyu"                        "Myeik"                          
[77] "Yangon (East)"                   "Yangon (North)"                 
[79] "Yangon (South)"                  "Yangon (West)"                  

There are 80 admin2 levels or districts in mmr_shp_mimu_2

Lets compare with our admin2 levels in our main dataset ACLED_MMR

unique_acled_regions2 <- unique(ACLED_MMR$admin2)

unique_acled_regions2

 [1] "Maungdaw"                        "Bago"                           
 [3] "Shwebo"                          "Kyaukme"                        
 [5] "Pyinoolwin"                      "Muse"                           
 [7] "Sittwe"                          "Yinmarbin"                      
 [9] "Thaton"                          "Yangon-North"                   
[11] "Pa-O Self-Administered Zone"     "Hpapun"                         
[13] "Kyaukpyu"                        "Yangon-West"                    
[15] "Mongmit"                         "Bhamo"                          
[17] "Mrauk-U"                         "Yangon-East"                    
[19] "Yangon-South"                    "Monywa"                         
[21] "Gangaw"                          "Pathein"                        
[23] "Katha"                           "Taungoo"                        
[25] "Kanbalu"                         "Lashio"                         
[27] "Mawlamyine"                      "Myitkyina"                      
[29] "Kawkareik"                       "Loilen"                         
[31] "Mandalay"                        "Kawlin"                         
[33] "Kyaukse"                         "Magway"                         
[35] "Meiktila"                        "Pakokku"                        
[37] "Taunggyi"                        "Tamu"                           
[39] "Nay Pyi Taw"                     "Mohnyin"                        
[41] "Kale"                            "Det Khi Na"                     
[43] "Myingyan"                        "Loikaw"                         
[45] "Matupi"                          "Pyay"                           
[47] "Sagaing"                         "Myeik"                          
[49] "Dawei"                           "Thayarwady"                     
[51] "Thandwe"                         "Mawlaik"                        
[53] "Bawlake"                         "Pyapon"                         
[55] "Hinthada"                        "Thayet"                         
[57] "Pa Laung Self-Administered Zone" "Mindat"                         
[59] "Hkamti"                          "Kokang Self-Administered Zone"  
[61] "Hpa-An"                          "Danu Self-Administered Zone"    
[63] "Myawaddy"                        "Maubin"                         
[65] "Hakha"                           "Falam"                          
[67] "Minbu"                           "Monghsat"                       
[69] "Puta-O"                          "Hopang"                         
[71] "Nyaung-U"                        "Kawthoung"                      
[73] "Yamethin"                        "Yangon"                         
[75] "Myaungmya"                       "Mong Pawk (Wa SAD)"             
[77] "Oke Ta Ra"                       "Matman"                         
[79] "Kengtung"                        "Naga Self-Administered Zone"    
[81] "Labutta"                         "Langkho"                        
[83] "Tachileik"                      

I will write a simple function below to identify the discrepancies between the shape file and our state/district names in our main dataset.

# Find the unique region names that are in 'unique_acled_regions2' but not in 'unique_regions_mimu2'

mismatched_admin2 <- setdiff(unique_acled_regions2, unique_regions_mimu2)

if (length(mismatched_admin2) > 0) {
  print("The following region names from 'acled_mmr' do not match any in 'mimu2':")
  print(mismatched_admin2)
} else {
  print("All unique region names in 'acled_mmr' match the unique region names in 'mimu2.'")
}

[1] "The following region names from 'acled_mmr' do not match any in 'mimu2':"
[1] "Yangon-North"       "Yangon-West"        "Yangon-East"       
[4] "Yangon-South"       "Nay Pyi Taw"        "Yangon"            
[7] "Mong Pawk (Wa SAD)"

Lets harmonize the names in both data files. I will re-save it to a new data set called ACLED_MMR_1

Fixing our admin 1 names.

ACLED_MMR_1 <- ACLED_MMR %>%
  mutate(admin1 = case_when(
    admin1 == "Bago-East" ~ "Bago (East)",
    admin1 == "Bago-West" ~ "Bago (West)",
    admin1 == "Shan-North" ~ "Shan (North)",
    admin1 == "Shan-South" ~ "Shan (South)",
    admin1 == "Shan-East" ~ "Shan (East)",
    TRUE ~ as.character(admin1)
  ))

Fixing our admin 2 names.

ACLED_MMR_1 <- ACLED_MMR_1 %>%
  mutate(admin2 = case_when(
    admin2 == "Yangon-East" ~ "Yangon (East)",
    admin2 == "Yangon-West" ~ "Yangon (West)",
    admin2 == "Yangon-North" ~ "Yangon (North)",
    admin2 == "Yangon-South" ~ "Yangon (South)",
    admin2 == "Mong Pawk (Wa SAD)" ~ "Tachileik",
    admin2 == "Nay Pyi Taw" ~ "Det Khi Na",
    admin2 == "Yangon" ~ "Yangon (West)",
    TRUE ~ as.character(admin2)
  ))

Checking if our changes are successful.

# Get unique admin 2 district names from 'ACLED_MMR_1'
unique_acled_regions2 <- unique(ACLED_MMR_1$admin2)

# Get unique district names from 'mmr_shp_mimu_2'
unique_map_regions_mimu2 <- unique(mmr_shp_mimu_2$DT)

# Find the unique district names that are in 'unique_acled_regions2' but not in 'unique_map_regions_mimu2'

mismatched_regions2 <- setdiff(unique_acled_regions2, unique_map_regions_mimu2)

if (length(mismatched_regions2) > 0) {
  print("The following district names from 'acled_mmr_1' do not match any in 'mmr_shp_mimu_2':")
  print(mismatched_regions2)
} else {
  print("All unique district names in 'acled_mmr_1' match the unique district names in 'mmmr_shp_mimu_2.'")
}

[1] "All unique district names in 'acled_mmr_1' match the unique district names in 'mmmr_shp_mimu_2.'"

Lets do a sample plot to see how our country map looks like at the admin2 (districts) level.

plot(mmr_shp_mimu_2)

4.3.2 Data Wrangling

For the purposes of plotting choropleth maps, I will first create attributes subsets to summarise the number of incidents and fatalities, grouped by year, admin region, event type and sub event type.

Data2 <- ACLED_MMR_1 %>%
    group_by(year, admin2, event_type, sub_event_type) %>%
    summarise(Incidents = n(),
              Fatalities = sum(fatalities, na.rm = TRUE)) %>%
              
    ungroup()

Checking the total sum of incidents and fatalities

total_incidents2 <- sum(Data2$Incidents)
total_fatalities2 <- sum(Data2$Fatalities)

total_incidents2

[1] 57198

total_fatalities2

[1] 57593

Next, I will do a spatial join between my shape files and attribute files

ACLED_MMR_admin2 <- left_join(mmr_shp_mimu_2, Data2,
                            by = c("DT" = "admin2"))

Removing the variables I don’t require.

ACLED_MMR_admin2 <- ACLED_MMR_admin2 %>%
                      select(-OBJECTID, -ST, -ST_PCODE)

class(ACLED_MMR_admin2)

[1] "sf"         "data.frame"

Next, I will just double check that total sum of incidents and fatalities in our SF files are correct as per our original datasets.

total_incidents_check <- sum(ACLED_MMR_admin2$Incidents)
total_fatalities_check <- sum(ACLED_MMR_admin2$Fatalities)

total_incidents_check 

[1] 57198

total_fatalities_check

[1] 57593

4.4 Choropleth Maps

For the module on Cluster & Outlier Analysis (LISA), for the first tab, I will be plotting both Choropleth and Proportional Symbol Maps, as part of the initial Interactive Exploratory Analysis..

4.4.1 Choropleth map of Incidents & Fatalities by Admin2 level (by District)

The below codes will be used to create the choropleth maps.

Fatalities in Battles in 2023, by Districts (Quantile)

tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +

ACLED_MMR_admin2 %>%
  filter(year == 2023, event_type == "Battles") %>%
  tm_shape() +
  tm_fill("Fatalities",
          n = 5,
          style = "quantile",
          palette = "Reds") +
  tm_borders(alpha = 0.5)

Incidents of Violence against civilians in 2021, by Districts (Equal)

tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +

ACLED_MMR_admin2 %>%
  filter(year == 2021, event_type == "Violence against civilians") %>%
  tm_shape() +
  tm_fill("Incidents",
          n = 5,
          style = "equal",
          palette = "Reds") +
  tm_borders(alpha = 0.5)

Adding Interactivity by using tmap_leaftlet()

data_filtered <- ACLED_MMR_admin2 %>%
  filter(year == 2022, event_type == "Battles")

tm_map <- tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +
  
  tm_shape(data_filtered) +
  tm_fill(col = "Incidents", n = 5, style = "equal", palette = "Reds", title = "Incidents") +
  tm_borders(alpha = 0.5)

tmap_leaflet(tm_map)

Parameters and Output to be exposed to Shiny

From the codes above, below are the variables we can expose as user inputs:-

• Year

• Event type (Battles, Violence against civilians, protests, riots, explosions/remote violence)

• Count type: number of Incidents or Fatalities

• Data classification type: eg quantile, equal, jenks, kmeans, pretty etc

4.5 Proportional Symbol Maps

Proportional symbol maps (also known as graduate symbol maps) are a class of maps that use the visual variable of size to represent differences in the magnitude of a discrete, abruptly changing phenomenon, e.g. counts of incidents, fatalities.

First I will convert the ACLED_MMR_1 data set to become an sf object called conflict_sf.

# Convert conflict data to an sf object
conflict_sf <- st_as_sf(ACLED_MMR_1, coords = c("longitude", "latitude"), crs = 4326)

class(conflict_sf)

[1] "sf"         "tbl_df"     "tbl"        "data.frame"

conflict_sf

Simple feature collection with 57198 features and 33 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 92.199 ymin: 9.9824 xmax: 100.3576 ymax: 27.731
Geodetic CRS:  WGS 84
# A tibble: 57,198 × 34
   event_id_cnty event_date        year time_precision disorder_type  event_type
 * <chr>         <chr>            <dbl>          <dbl> <chr>          <chr>     
 1 MMR58558      16 February 2024  2024              1 Political vio… Battles   
 2 MMR58559      16 February 2024  2024              1 Political vio… Battles   
 3 MMR58443      16 February 2024  2024              2 Political vio… Violence …
 4 MMR58502      16 February 2024  2024              1 Political vio… Violence …
 5 MMR58507      16 February 2024  2024              1 Political vio… Explosion…
 6 MMR58508      16 February 2024  2024              1 Political vio… Explosion…
 7 MMR58547      16 February 2024  2024              1 Strategic dev… Strategic…
 8 MMR58560      16 February 2024  2024              1 Political vio… Battles   
 9 MMR58589      16 February 2024  2024              2 Political vio… Violence …
10 MMR58648      16 February 2024  2024              1 Political vio… Explosion…
# ℹ 57,188 more rows
# ℹ 28 more variables: sub_event_type <chr>, actor1 <chr>, assoc_actor_1 <chr>,
#   inter1 <dbl>, actor2 <chr>, assoc_actor_2 <chr>, inter2 <dbl>,
#   interaction <dbl>, civilian_targeting <chr>, iso <dbl>, region <chr>,
#   country <chr>, admin1 <chr>, admin2 <chr>, admin3 <chr>, location <chr>,
#   geo_precision <dbl>, source <chr>, source_scale <chr>, notes <chr>,
#   fatalities <dbl>, tags <chr>, timestamp <dbl>, population_1km <dbl>, …

Next, I create subsets for each event type, each subset will inherit the SF object characteristic.

Battles <- filter(conflict_sf, event_type == "Battles")
Violence_CV <- filter(conflict_sf, event_type == "Violence against civilians")
Protests <- filter(conflict_sf, event_type == "Protests")
Riots <- filter(conflict_sf, event_type == "Riots")
Explosions <- filter(conflict_sf, event_type == "Explosions/Remote violence")
Strategic_dev <- filter(conflict_sf, event_type == "Strategic developments")

class(Battles)

[1] "sf"         "tbl_df"     "tbl"        "data.frame"

4.5.1 Visualising the location of conflict events

Using the leaflet package, I will use the Geometry points from our SF data sets to plot the points of event types in the maps.

In this case, I will use admin 1 level regions, to achieve a better aesthetics for users.

This is because, visually dividing the country map into more smaller districts (admin2) would likely make the map look “too busy”.

Battles from 2010 to present

scaleFactor <- 2  

leaflet(Battles) %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  addPolygons(data = mmr_shp_mimu_1, color = "#444444", weight = 1, fillOpacity = 0.5) %>% # Adding borders
  
  addCircleMarkers(popup = ~paste("Event: Battles<br>State/Region:", admin1, 
                                  "<br>Actor1:", actor1, "<br>Actor2:", actor2,
                                  "<br>Year:", year, "<br>Fatalities:", fatalities),
                   radius = ~sqrt(fatalities) * scaleFactor,
                   fillColor = "red", fillOpacity = 0.4, color = "#FFFFFF", weight = 1) %>% 
  setView(lng = 96.1603, lat = 19.745, zoom = 6)

Violence against civilians from 2010 to present

scaleFactor <- 2  

leaflet(Violence_CV) %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  addPolygons(data = mmr_shp_mimu_1, color = "#444444", weight = 1, fillOpacity = 0.5) %>% # Adding borders
  
  addCircleMarkers(popup = ~paste("Event: Violence on Civillians<br>State/Region:", admin1, 
                                  "<br>Actor1:", actor1, "<br>Actor2:", actor2,
                                  "<br>Year:", year, "<br>Fatalities:", fatalities),
                   radius = ~sqrt(fatalities) * scaleFactor,
                   fillColor = "red", fillOpacity = 0.4, color = "#FFFFFF", weight = 1) %>% 
  setView(lng = 96.1603, lat = 19.745, zoom = 6)

Protests from 2010 to present

scaleFactor <- 2  

leaflet(Protests) %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  addPolygons(data = mmr_shp_mimu_1, color = "#444444", weight = 1, fillOpacity = 0.5) %>% # Adding borders
  
  addCircleMarkers(popup = ~paste("Event: Protests<br>State/Region:", admin1, 
                                  "<br>Actor1:", actor1, "<br>Actor2:", actor2,
                                  "<br>Year:", year, "<br>Fatalities:", fatalities),
                   radius = ~sqrt(fatalities) * scaleFactor,
fillColor = "red", fillOpacity = 0.4, color = "#FFFFFF", weight = 1) %>% 
  setView(lng = 96.1603, lat = 19.745, zoom = 6)

Parameters and Output to be exposed to Shiny

The plots above are sufficiently interactive as users can hover to any “circle” to get more information on the event and location.

In terms of enabling user inputs, I will expose it to accept user input: Year and Event type.

4.6 UI design - Part 1

As mentioned, both plots above will be used to populate our first tab in our app’s Cluster & Outlier Analysis module.

The motivation here is to provide users with a Visual Introduction of the state of Armed Conflict events in the country (Myanmar), before they proceed to other tabs that has more emphasis on spacial autocorrelation and cluster analysis.

Below is a visual of the prototype for this page.

Functionality and interactivity

• The Proportional symbol map will allow users to select the specific year and event type. Users will be able to see where conflicts events have happened in the country across different years.

• The Choropleth map will be able to visualise the “intensity” of the events through use of color gradient. Users will be able to select the specific year, event type, count type (incidents or fatalities) and the data classification method.

Additional considerations

For the symbol map, I have populated the tool tip to communicate the event type, year, region name, number of fatalities and the actors involved in the event.

Instead of implementing a global filter which can enable users to filter and affect both plots, I have chosen to apply seperate filter selections for each plot. This is to enable users to explore different data independently. For example, users can explore “Battles” in 2021 via the symbol map and “Violence against civilians” in 2021 via the choropleth map.

4.7 Data Preparation for Spatial Analysis

For the next part of our UI, we will be exploring Spatial statistics.

Specifically I will be deriving and visualising the Local Moran’s I statistics.

First, I will create subsets of our Events happening in admin region 1 & 2, summarized with the number(count) of incidents and fatalities.

Events1 <- ACLED_MMR_1 %>%
    group_by(year, admin1, event_type) %>%
    summarise(Incidents = n(),
              Fatalities = sum(fatalities, na.rm = TRUE)) %>%
              
    ungroup()


Events2 <- ACLED_MMR_1 %>%
    group_by(year, admin2, event_type) %>%
    summarise(Incidents = n(),
              Fatalities = sum(fatalities, na.rm = TRUE)) %>%
              
    ungroup()

Next, I will perform a relational join to update our admin 1 and admin 2 level shape files with attributes fields of the above event related data.

Events_admin1 <- left_join(mmr_shp_mimu_1, Events1,
                            by = c("ST" = "admin1"))

Events_admin2 <- left_join(mmr_shp_mimu_2, Events2,
                            by = c("DT" = "admin2"))

Removing the variables I don’t need.

Events_admin2 <- Events_admin2 %>%
                      select(-OBJECTID, -ST, -ST_PCODE)

class(Events_admin2)

[1] "sf"         "data.frame"

st_geometry(Events_admin2)

Geometry set for 2748 features 
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS:  WGS 84
First 5 geometries:

Next, I will create a subset of the event type and year.

For the puposes of this exercise, which is to test the code outputs, I will create a subset to analyse for Event type = Battles, in the year 2023.

I will name the object as Battles_2023, eventually this object file will be coded and named generically, in our app and be used to “carry” users selection (eg event types and year).

4.7.1 Filtering the Event and Year (Event type = Battles, in 2023)

The below subset will serve as our reference data subset for our subsequent codes.

Battles_2023 <- Events_admin2 %>%
  filter(year == 2023, event_type == "Battles")

4.7.2 Computing Contiguity Spatial Weights

Before we can compute any spatial statistics, we need to construct spatial weights of the study area.

The spatial weights is used to define the neighbourhood relationships between the geographical units (i.e. admin2) in the study area (Myanmar).

In the code below, poly2nb() of spdep package is used to compute contiguity weight matrices for the study area. This function builds a neighbours list based on regions with contiguous boundaries.

By default this function will return a list of first order neighbours using the Queen criteria.

However, we can also pass a “queen” argument that takes TRUE or FALSE as options.

wm_q <- poly2nb(Battles_2023, 
                queen=TRUE)
summary(wm_q)

Neighbour list object:
Number of regions: 74 
Number of nonzero links: 356 
Percentage nonzero weights: 6.501096 
Average number of links: 4.810811 
Link number distribution:

 1  2  3  4  5  6  7  8 10 
 3  8 11 11 15  9  9  6  2 
3 least connected regions:
2 16 58 with 1 link
2 most connected regions:
19 56 with 10 links

The summary report above shows that there are 74 area units for this subset (Battles occurring in 2023).

There are 2 most connected area units with 10 neighbours, and there are 3 area units with only 1 neighbour.

4.7.3 Row-standardised weights matrix

Next, we assign weights to each neighboring polygon. In our case, each neighboring polygon will be assigned equal weight (style=“W”). This is accomplished by assigning the fraction 1/(#ofneighbors) to each neighboring admin2 (district) and then summing the weighted income values.

This has one drawback in that polygons along the edges of the study area will base their lagged values on fewer polygons and thus potentially over or under estimating the true nature of the spatial autocorrelation in the data.

However, for this example, I will stick with the style=“W” option for simplicity’s sake. Other more robust options are available, notably style=“B”.

rswm_q <- nb2listw(wm_q, 
                   style="W", 
                   zero.policy = TRUE)
rswm_q

Characteristics of weights list object:
Neighbour list object:
Number of regions: 74 
Number of nonzero links: 356 
Percentage nonzero weights: 6.501096 
Average number of links: 4.810811 

Weights style: W 
Weights constants summary:
   n   nn S0       S1       S2
W 74 5476 74 36.81702 310.0455

4.8 Cluster and Outlier Analysis

While GLOBAL Moran’s I score and the Geary’s C ratio can tell us whether specific event types (Battles, Explosions, Protests, Riots, Violence against civilians) tends to cluster or not on the map, it does not provide any information on the distribution of spatial dependence of Events types, and is unable to identify the location of hotspots and clusters.

For that, we require the use of more localized methods - Anselin’s Moran Scatterplot and the Local Indicator of Spatial Autocorrelation (LISA) method.

Local Indicators of Spatial Association or LISA are statistics that evaluate the existence of clusters in the spatial arrangement of a given variable.

For example, in this analysis, we are studying if there are areas that have higher or lower incidents of a specific Event type (Battles) than is to be expected by chance alone, ie the values occurring are above or below those of a random distribution in space.

4.8.1 Computing Local Moran’s I

To compute local Moran’s I, the localmoran() function of spdep will be used. It computes Ii values, given a set of zi values and a list object providing neighbour weighting information for the polygon associated with the zi values.

fips <- order(Battles_2023$DT)
localMI <- localmoran(Battles_2023$Incidents, rswm_q)
head(localMI)

           Ii          E.Ii       Var.Ii       Z.Ii Pr(z != E(Ii))
1  0.46274887 -9.006432e-03 0.1247560902  1.3356292      0.1816705
2  0.71317847 -1.016877e-02 0.7448368248  0.8381394      0.4019524
3  0.69218038 -9.386041e-03 0.3392457987  1.2045132      0.2283913
4  0.68518101 -9.386041e-03 0.3392457987  1.1924960      0.2330668
5  0.00627746 -1.483579e-05 0.0001702656  0.4822206      0.6296493
6 -0.23004674 -4.814215e-03 0.0464274459 -1.0453066      0.2958813

localmoran() function returns a matrix of values whose columns are:

• Ii: the local Moran’s I statistics

• E.Ii: the expectation of local moran statistic under the randomisation hypothesis

• Var.Ii: the variance of local moran statistic under the randomisation hypothesis

• Z.Ii:the standard deviate of local moran statistic

• Pr(): the p-value of local moran statistic

The code below lists the content of the local Moran matrix derived by using printCoefmat().

printCoefmat(data.frame(
  localMI[fips,], 
  row.names=Battles_2023$DT[fips]),
  check.names=FALSE)

                                         Ii        E.Ii      Var.Ii        Z.Ii
Bago                             6.2775e-03 -1.4836e-05  1.7027e-04  4.8222e-01
Bawlake                         -1.7909e-02 -4.6941e-04  1.1252e-02 -1.6441e-01
Bhamo                            2.4621e-01 -1.7367e-03  1.9897e-02  1.7577e+00
Danu Self-Administered Zone      3.3657e-01 -8.2707e-03  1.9670e-01  7.7752e-01
Dawei                            9.4559e-01 -1.4130e-02  5.0825e-01  1.3462e+00
Det Khi Na                       2.3024e-01 -6.5686e-03  6.3234e-02  9.4170e-01
Falam                           -7.5203e-03 -2.4381e-03  5.8326e-02 -2.1044e-02
Gangaw                           6.2609e-01 -6.9289e-03  6.6679e-02  2.4514e+00
Hakha                           -1.1347e-02 -6.1003e-05  1.0815e-03 -3.4318e-01
Hinthada                         4.6275e-01 -9.0064e-03  1.2476e-01  1.3356e+00
Hkamti                          -5.3716e-02 -3.0598e-03  3.5009e-02 -2.7074e-01
Hopang                          -2.3647e-01 -9.7735e-03  3.5311e-01 -3.8150e-01
Hpa-An                          -1.3231e-01 -3.0598e-03  1.9751e-02 -9.1968e-01
Hpapun                          -4.1969e-02 -4.2526e-03  5.9189e-02 -1.5503e-01
Kale                             1.2532e-01 -7.7384e-04  6.4571e-03  1.5692e+00
Kanbalu                         -3.2905e-02 -3.4002e-05  3.9022e-04 -1.6640e+00
Katha                           -8.8073e-03 -1.8682e-02  1.5309e-01  2.5238e-02
Kawkareik                        2.9062e-02 -1.3663e-02  3.2318e-01  7.5156e-02
Kawlin                          -7.2541e-02 -1.4063e-03  3.3678e-02 -3.8762e-01
Kawthoung                       -6.7102e-01 -3.2826e-03  2.4212e-01 -1.3570e+00
Kokang Self-Administered Zone    9.6449e-04 -1.1447e-08  2.7452e-07  1.8408e+00
Kyaukme                         -6.9939e-02 -9.0471e-03  8.6877e-02 -2.0659e-01
Kyaukpyu                         4.7922e-01 -6.8933e-03  1.2137e-01  1.3953e+00
Kyaukse                         -1.4472e-01 -9.0064e-03  7.4533e-02 -4.9713e-01
Labutta                          7.1318e-01 -1.0169e-02  7.4484e-01  8.3814e-01
Langkho                         -9.6769e-03 -8.6347e-03  2.0528e-01 -2.3004e-03
Lashio                           2.0219e-01 -4.2805e-03  4.8917e-02  9.3352e-01
Loikaw                          -4.7652e-01 -2.3869e-02  3.2568e-01 -7.9318e-01
Loilen                          -2.8804e-02 -5.6413e-03  7.8408e-02 -8.2718e-02
Magway                           9.4077e-02 -5.9426e-03  4.9330e-02  4.5033e-01
Mandalay                        -2.2955e-01 -3.0598e-03  7.3153e-02 -8.3742e-01
Matupi                          -1.7986e-02 -3.2117e-04  4.4878e-03 -2.6370e-01
Maungdaw                         1.2575e-01 -2.6375e-03  6.3084e-02  5.1117e-01
Mawlaik                         -3.2557e-02 -2.6375e-03  3.0190e-02 -1.7220e-01
Mawlamyine                       2.4766e-01 -3.3072e-03  5.8440e-02  1.0381e+00
Meiktila                         3.8976e-01 -8.2707e-03  7.9484e-02  1.4118e+00
Minbu                           -1.9321e-01 -6.2516e-03  6.0203e-02 -7.6195e-01
Mindat                           2.2211e-02 -8.7518e-04  1.5503e-02  1.8541e-01
Mohnyin                          1.1002e-01 -8.8788e-04  1.5727e-02  8.8440e-01
Mongmit                         -4.9751e-01 -8.6347e-03  1.1965e-01 -1.4133e+00
Monywa                           3.7540e+00 -3.3925e-02  5.8106e-01  4.9692e+00
Mrauk-U                          2.1574e-01 -3.9983e-03  4.5705e-02  1.0278e+00
Muse                             2.6578e-01 -1.6313e-01  2.4204e+00  2.7569e-01
Myawaddy                         1.8035e-02 -6.4391e-05  2.3492e-03  3.7341e-01
Myeik                            3.6057e-01 -2.5739e-02  9.1496e-01  4.0386e-01
Myingyan                         4.3732e-02 -1.0008e-04  1.3987e-03  1.1720e+00
Myitkyina                       -1.9413e-01 -6.2855e-03  8.7306e-02 -6.3572e-01
Naga Self-Administered Zone     -8.8212e-02 -1.0169e-02  3.6725e-01 -1.2878e-01
Nyaung-U                        -5.4115e-01 -8.6347e-03  1.5176e-01 -1.3669e+00
Oke Ta Ra                        5.7519e-01 -9.0064e-03  1.5824e-01  1.4686e+00
Pa-O Self-Administered Zone      1.9941e-01 -5.6413e-03  5.4359e-02  8.7950e-01
Pa Laung Self-Administered Zone -5.3309e-01 -5.0623e-03  7.0402e-02 -1.9901e+00
Pakokku                          1.9341e+00 -2.2766e-01  1.4683e+00  1.7840e+00
Pathein                          6.9218e-01 -9.3860e-03  3.3925e-01  1.2045e+00
Puta-O                          -4.8052e-01 -6.8933e-03  5.0659e-01 -6.6543e-01
Pyapon                           6.8518e-01 -9.3860e-03  3.3925e-01  1.1925e+00
Pyay                             1.0545e-01 -1.2618e-03  1.7615e-02  8.0407e-01
Pyinoolwin                       4.6151e-01 -2.7025e-02  2.1958e-01  1.0426e+00
Sagaing                          9.5824e-01 -1.0212e-02  9.7948e-02  3.0944e+00
Shwebo                           2.1412e+00 -5.0075e-02  5.4593e-01  2.9657e+00
Sittwe                           2.3135e-01 -3.0598e-03  1.1130e-01  7.0266e-01
Tamu                             3.2174e-02 -1.8902e-04  3.3506e-03  5.5911e-01
Taunggyi                        -1.0168e-01 -1.1395e-03  7.3697e-03 -1.1711e+00
Taungoo                         -2.3005e-01 -4.8142e-03  4.6427e-02 -1.0453e+00
Thandwe                          5.6456e-01 -1.0169e-02  1.4069e-01  1.5322e+00
Thaton                          -6.7767e-02 -3.0835e-03  5.4499e-02 -2.7708e-01
Thayarwady                       9.0973e-03 -1.4836e-05  2.0737e-04  6.3277e-01
Thayet                           2.9530e-01 -5.3479e-03  5.1546e-02  1.3242e+00
Yamethin                         4.3920e-01 -9.7735e-03  1.3528e-01  1.2207e+00
Yangon (East)                    6.1100e-01 -7.5664e-03  1.8008e-01  1.4576e+00
Yangon (North)                   4.4954e-01 -9.3860e-03  1.0671e-01  1.4049e+00
Yangon (South)                   5.2023e-01 -8.6347e-03  1.1965e-01  1.5289e+00
Yangon (West)                    6.6585e-01 -9.7735e-03  2.3209e-01  1.4024e+00
Yinmarbin                        4.5788e+00 -9.9115e-02  1.2481e+00  4.1872e+00
                                Pr.z....E.Ii..
Bago                                    0.6296
Bawlake                                 0.8694
Bhamo                                   0.0788
Danu Self-Administered Zone             0.4369
Dawei                                   0.1782
Det Khi Na                              0.3463
Falam                                   0.9832
Gangaw                                  0.0142
Hakha                                   0.7315
Hinthada                                0.1817
Hkamti                                  0.7866
Hopang                                  0.7028
Hpa-An                                  0.3577
Hpapun                                  0.8768
Kale                                    0.1166
Kanbalu                                 0.0961
Katha                                   0.9799
Kawkareik                               0.9401
Kawlin                                  0.6983
Kawthoung                               0.1748
Kokang Self-Administered Zone           0.0656
Kyaukme                                 0.8363
Kyaukpyu                                0.1629
Kyaukse                                 0.6191
Labutta                                 0.4020
Langkho                                 0.9982
Lashio                                  0.3506
Loikaw                                  0.4277
Loilen                                  0.9341
Magway                                  0.6525
Mandalay                                0.4024
Matupi                                  0.7920
Maungdaw                                0.6092
Mawlaik                                 0.8633
Mawlamyine                              0.2992
Meiktila                                0.1580
Minbu                                   0.4461
Mindat                                  0.8529
Mohnyin                                 0.3765
Mongmit                                 0.1576
Monywa                                  0.0000
Mrauk-U                                 0.3040
Muse                                    0.7828
Myawaddy                                0.7088
Myeik                                   0.6863
Myingyan                                0.2412
Myitkyina                               0.5250
Naga Self-Administered Zone             0.8975
Nyaung-U                                0.1716
Oke Ta Ra                               0.1419
Pa-O Self-Administered Zone             0.3791
Pa Laung Self-Administered Zone         0.0466
Pakokku                                 0.0744
Pathein                                 0.2284
Puta-O                                  0.5058
Pyapon                                  0.2331
Pyay                                    0.4214
Pyinoolwin                              0.2972
Sagaing                                 0.0020
Shwebo                                  0.0030
Sittwe                                  0.4823
Tamu                                    0.5761
Taunggyi                                0.2415
Taungoo                                 0.2959
Thandwe                                 0.1255
Thaton                                  0.7817
Thayarwady                              0.5269
Thayet                                  0.1854
Yamethin                                0.2222
Yangon (East)                           0.1449
Yangon (North)                          0.1601
Yangon (South)                          0.1263
Yangon (West)                           0.1608
Yinmarbin                               0.0000

4.8.2 Mapping the Local Moran’s I

Before mapping the local Moran’s I map, I will need to append the local Moran’s I dataframe (i.e. localMI) onto the Battles_2023’s SF DataFrame.

Battles_2023.localMI <- cbind(Battles_2023,localMI) %>%
  rename(Pr.Ii = Pr.z....E.Ii..)

Battles_2023.localMI

Simple feature collection with 74 features and 13 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 99.66532 ymax: 28.54554
Geodetic CRS:  WGS 84
First 10 features:
           DT   DT_PCODE    DT_MMR PCode_V year event_type Incidents Fatalities
1    Hinthada MMR017D002  ဟင်္သာတခရိုင်     9.4 2023    Battles         4          3
2     Labutta MMR017D004  လပွတ္တာခရိုင်     9.4 2023    Battles         1          1
3     Pathein MMR017D001    ပုသိမ်ခရိုင်     9.4 2023    Battles         3          1
4      Pyapon MMR017D006   ဖျာပုံခရိုင်     9.4 2023    Battles         3          2
5        Bago MMR007D001    ပဲခူးခရိုင်     9.4 2023    Battles        50        270
6     Taungoo MMR007D002  တောင်ငူခရိုင်     9.4 2023    Battles        87        493
7        Pyay MMR008D001    ပြည်ခရိုင်     9.4 2023    Battles        34         59
8  Thayarwady MMR008D002 သာယာဝတီခရိုင်     9.4 2023    Battles        50         93
9       Falam MMR004D001   ဖလမ်းခရိုင်     9.4 2023    Battles        27         89
10      Hakha MMR004D003 ဟားခါးခရိုင်     9.4 2023    Battles        48        210
             Ii          E.Ii       Var.Ii        Z.Ii     Pr.Ii
1   0.462748867 -9.006432e-03 0.1247560902  1.33562921 0.1816705
2   0.713178468 -1.016877e-02 0.7448368248  0.83813940 0.4019524
3   0.692180375 -9.386041e-03 0.3392457987  1.20451317 0.2283913
4   0.685181011 -9.386041e-03 0.3392457987  1.19249602 0.2330668
5   0.006277460 -1.483579e-05 0.0001702656  0.48222056 0.6296493
6  -0.230046740 -4.814215e-03 0.0464274459 -1.04530664 0.2958813
7   0.105454317 -1.261774e-03 0.0176145487  0.80407054 0.4213562
8   0.009097304 -1.483579e-05 0.0002073682  0.63277490 0.5268807
9  -0.007520292 -2.438086e-03 0.0583263538 -0.02104359 0.9832109
10 -0.011346560 -6.100301e-05 0.0010814666 -0.34317564 0.7314663
                         geometry
1  MULTIPOLYGON (((95.12637 18...
2  MULTIPOLYGON (((95.04462 15...
3  MULTIPOLYGON (((94.27572 15...
4  MULTIPOLYGON (((95.20798 15...
5  MULTIPOLYGON (((95.90674 18...
6  MULTIPOLYGON (((96.17964 19...
7  MULTIPOLYGON (((95.70458 19...
8  MULTIPOLYGON (((95.85173 18...
9  MULTIPOLYGON (((93.36931 24...
10 MULTIPOLYGON (((93.35213 23...

4.8.3 Mapping Local Moran’s I values

Using choropleth mapping functions of tmap package, we can plot the local Moran’s I values by using the code below.

tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +

tm_shape(Battles_2023.localMI) +
  tm_fill(col = "Ii", 
          style = "pretty",
          palette = "RdBu",
          title = "local moran statistics") +
  tm_borders(alpha = 0.5)

4.8.4 Mapping Local Moran’s I p-values

The choropleth above shows there is evidence for both positive and negative Ii values. However, we will also need to consider the p-values for each of these values.

The code below produces a choropleth map of Moran’s I p-values by using functions of tmap package.

tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +

tm_shape(Battles_2023.localMI) +
  tm_fill(col = "Pr.Ii", 
          breaks=c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
          palette="-Blues", 
          title = "local Moran's I p-values") +
  tm_borders(alpha = 0.5)

Parameters and Output to be exposed to Shiny

In terms of enabling user inputs, I will expose it to accept user input: Year and Event type.

4.8.6 Data table for the Moran’s I values

For the sake of readability, it may also be a good idea to add a data table of the values, for users to make sense of both maps.

The below code will be used to generate the data table.

datatable(Battles_2023.localMI)

4.9 UI design - Part 2

As mentioned, both plots and the data table above will be used to populate our second tab in our app’s Cluster and Outlier Analysis module.

The motivation here is to provide users with a visualisation of the Local Moran’s I statistics before they proceed to the next tab that has the Moran’s scatter plot and Lisa Cluster map.

The below is a visual of the prototype.

Functionality and interactivity

• Users will only need to select the year and the event type once. Both maps, along with the data table will be updated upon user selection.

Additional considerations

Here, I have chosen to implement a single point for users to filter and select, as users are likely to want to see all 3 plots communicating the same statistics.

4.10 Creating the LISA cluster map

For the next tab of our Cluster and Outlier Analysis module, I will be creating the Moran Scatter plot and the LISA cluster map.

The LISA Cluster Map shows the significant locations color coded by type of spatial autocorrelation.

4.10.1 Plotting the Moran Scatter plot

The Moran scatterplot is an illustration of the relationship between the values of the chosen attribute at each location and the average value of the same attribute at neighboring locations.

The code below plots the Moran scatterplot of Battles in 2023 by using moran.plot() of spdep.

nci <- moran.plot(Battles_2023$Incidents, rswm_q,
                  labels=as.character(Battles_2023$DT), 
                  xlab="Battles_2023", 
                  ylab="Spatially Lagged Events,Year")

The plot is split in 4 quadrants. The top right corner belongs to areas that have high incidents of events and are surrounded by other areas that have higher than the average level/number of battles This is the high-high locations.

Note

The Moran scatterplot is divided into four areas, with each quadrant corresponding with one of four categories: (1) High-High (HH) in the top-right quadrant; (2) High-Low (HL) in the bottom right quadrant; (3) Low-High (LH) in the top-left quadrant; (4) Low- Low (LL) in the bottom left quadrant.

4.10.2 Plotting Moran scatterplot with standardised variable

First, I will use scale() to centre and scale the variable. Here centering is done by subtracting the mean (omitting NAs) the corresponding columns, and scaling is done by dividing the (centred) variable by their standard deviations.

Battles_2023$Z.Incidents <- scale(Battles_2023$Incidents) %>% 
  as.vector 

The as.vector() is added to the end is to make sure that the data type is a vector, that maps neatly into our dataframe.

Next, we plot the Moran scatterplot again by using the code below.

nci2 <- moran.plot(Battles_2023$Z.Incidents, rswm_q,
                   labels=as.character(Battles_2023$DT),
                   xlab="z-Battles in 2023", 
                   ylab="Spatially Lag z-Battles in 2023")

Note

1) High-High (HH): indicates high spatial correlation where incidents of Battles are clustered closely together.

2) High-Low (HL): where areas of high frequency of incidents of Battles occurred are located next to areas where there is low frequency of incidents of Battles occurred.

3) Low-High (LH): these are areas of low frequency of incidents where Battles occurred that are located next to areas where high frequency of Battles.

4) Low-Low (LL): these are clusters of low frequency of incidents of Battles occurred.

4.10.3 Preparing LISA map classes

According to Anselin (1995), LISA can be used to locate “hot spots” or local spatial clusters where the occurrence of Event types is statistically significant.

In addition to the four categories described in the Moran Scatterplot, the LISA analysis includes an additional category: (5) Insignificant: where there are no spatial autocorrelation or clusters where event types have occurred.

The code below shows the steps to prepare a LISA cluster map.

quadrant <- vector(mode="numeric",length=nrow(localMI))

Next, we derive the spatially lagged variable of interest (i.e. Incidents) and centre the spatially lagged variable around its mean.

Battles_2023$lag_Incidents <- lag.listw(rswm_q, Battles_2023$Incidents)
DV <- Battles_2023$lag_Incidents - mean(Battles_2023$lag_Incidents)     

This is followed by centering the local Moran’s around the mean.

LM_I <- localMI[,1] - mean(localMI[,1])    

Next, we will set a statistical significance level for the local Moran.

signif <- 0.05       

The code below define the low-low (1), low-high (2), high-low (3) and high-high (4) categories.

quadrant[DV <0 & LM_I>0] <- 1
quadrant[DV >0 & LM_I<0] <- 2
quadrant[DV <0 & LM_I<0] <- 3  
quadrant[DV >0 & LM_I>0] <- 4      

Lastly, we place non-significant Moran in the category 0.

quadrant[localMI[,5]>signif] <- 0

4.10.4 Plotting LISA Map

The below code is used to create the LISA map.

Battles_2023.localMI$quadrant <- quadrant
colors <- c("lightyellow", "#2c7bb6", "#abd9e9", "#fdae61", "#d7191c")
clusters <- c("insignificant", "low-low", "low-high", "high-low", "high-high")

tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +

tm_shape(Battles_2023.localMI) +
  tm_fill(col = "quadrant", 
          style = "cat", 
          palette = colors[c(sort(unique(quadrant)))+1], 
          labels = clusters[c(sort(unique(quadrant)))+1],
          popup.vars = c("")) +
  tm_view(set.zoom.limits = c(11,17)) +
  tm_borders(alpha=0.5)

Parameters and Output to be exposed to Shiny

In terms of enabling user inputs, I will expose it to accept user input: Year and Event type.

4.11 UI design - Part 3

As mentioned, both plots above will be used to populate our third tab in our app’s Cluster and Outlier Analysis module.

The motivation here is to provide users with a visualisation of statistically significant clusters of event types in order to help users understand how clusters may have developed and/or changed over time.

The below is a visual of the prototype.

Functionality and interactivity

• Both plots will enable users to select the specific year and event type.

• The Moran scatterplot enables users to see the statistically significant region names, by means of the quadrant and the annotations.

• The Lisa Cluster map adds value by showing the specific parts of the country which are categorized as statistically significant or insignificant.

Additional considerations

Instead of implementing a global filter which can enable users to filter and affect both plots, I have chosen to apply seperate filter selections for each plot. This is to enable users to explore different data independently. For example, users can explore “Battles” in 2021 via the scatter plot and “Violence against civilians” in 2021 via the LISA cluster map.

4.12 Hot & Cold Spot Area Analysis

In my last module, I will create plots to analyse the Hot and Cold spots across the country.

Beside detecting for clusters and outliers, Localised spatial statistics can be also used to detect hot spot and/or cold spot areas.

According to Lepers et al 2005, Aben et al 2012 and Isobe et al 2015; the term ‘hot spot’ has been used generically across disciplines to describe a region or value that is higher relative to its surroundings.

Unlike the previous section utilizing the Local Moran’s I statistics, here we will be utilizing the Getis and Ord’s G statistics (Getis and Ord, 1972; Ord and Getis, 1995).

This is an alternative spatial statistics for detecting spatial anomalies. It looks at neighbours within a defined proximity (distance) to identify where either high or low values clusters spatially.

Statistically significant hot-spots are recognised as areas of high values where other areas within a neighbourhood range also share high values too.

The workflow consists of 3 steps:

1. Deriving spatial weight matrix

2. Computing Gi statistics

3. Mapping Gi statistics

4.12.1 Deriving distance-based weight matrix

Whist the spatial autocorrelation in the previous section considered units which shared borders, for Getis-Ord we are defining neighbours based on distance.

There are two type of distance-based proximity matrix, they are:

• fixed distance weight matrix; and

• adaptive distance weight matrix.

4.12.1.1 Deriving the centroid

We will need points to associate with each polygon before we can make our connectivity graph.

We need the coordinates in a separate data frame for this to work. To do this we will use a mapping function. The mapping function applies a given function to each element of a vector and returns a vector of the same length.

Our input vector will be the geometry column of Battles_2023 dataset. Our function will be st_centroid(). We will be using map_dbl variation of map from the purrr package.

To get our longitude values we map the st_centroid() function over the geometry column of our Battles_2023 dataset and access the longitude value through double bracket notation [[]] and 1.

This allows us to get only the longitude, which is the first value in each centroid.

longitude <- map_dbl(Battles_2023$geometry, ~st_centroid(.x)[[1]])

class(longitude)

[1] "numeric"

longitude

 [1] 95.19035 94.99369 94.74008 95.53705 96.58767 96.34709 95.30186 95.75119
 [9] 93.71488 93.56355 93.22544 93.81953 97.16218 96.49805 97.41891 97.78205
[17] 97.42880 97.33311 97.45091 97.36866 98.22657 98.51991 94.18202 95.31595
[25] 94.45732 94.73063 95.05170 96.22693 92.48276 93.35447 92.90601 94.50024
[33] 95.47237 94.40944 95.49546 96.05576 95.53962 94.74604 95.26376 95.66196
[41] 95.64449 95.40740 94.38739 94.71565 98.96838 98.65407 97.13646 98.19561
[49] 96.66042 98.01229 97.21860 96.52970 98.02872 97.92760 97.03704 96.92822
[57] 98.47076 98.77515 96.16459 95.97341 95.54920 95.15564 96.25009 96.06338
[65] 97.84218 97.26005 96.16788 96.13921 93.92534 98.93770 96.24776 96.03014
[73] 96.28511 96.13431

We do the same for latitude by accessing the second value per centroid with [[2]].

latitude <- map_dbl(Battles_2023$geometry, ~st_centroid(.x)[[2]])

Now that we have latitude and longitude, we use cbind to put longitude and latitude into the same object.

coords <- cbind(longitude, latitude)

class(coords)

[1] "matrix" "array" 

coords

      longitude latitude
 [1,]  95.19035 17.87515
 [2,]  94.99369 16.16236
 [3,]  94.74008 16.86214
 [4,]  95.53705 16.17363
 [5,]  96.58767 17.74287
 [6,]  96.34709 18.81260
 [7,]  95.30186 18.76076
 [8,]  95.75119 18.10024
 [9,]  93.71488 23.42199
[10,]  93.56355 22.51244
[11,]  93.22544 21.57539
[12,]  93.81953 21.26333
[13,]  97.16218 24.25814
[14,]  96.49805 25.30215
[15,]  97.41891 26.01108
[16,]  97.78205 27.27547
[17,]  97.42880 18.91417
[18,]  97.33311 19.48546
[19,]  97.45091 17.75720
[20,]  97.36866 18.07659
[21,]  98.22657 16.00305
[22,]  98.51991 16.54942
[23,]  94.18202 21.80335
[24,]  95.31595 20.26548
[25,]  94.45732 20.38303
[26,]  94.73063 21.43135
[27,]  95.05170 19.43717
[28,]  96.22693 21.61640
[29,]  92.48276 21.02760
[30,]  93.35447 20.47696
[31,]  92.90601 20.40055
[32,]  94.50024 18.53159
[33,]  95.47237 25.34589
[34,]  94.40944 23.07978
[35,]  95.49546 23.32423
[36,]  96.05576 24.25502
[37,]  95.53962 23.93661
[38,]  94.74604 23.96919
[39,]  95.26376 22.25730
[40,]  95.66196 26.41610
[41,]  95.64449 21.98775
[42,]  95.40740 22.74644
[43,]  94.38739 24.17969
[44,]  94.71565 22.24865
[45,]  98.96838 23.02589
[46,]  98.65407 23.83364
[47,]  97.13646 22.50271
[48,]  98.19561 22.76825
[49,]  96.66042 23.45728
[50,]  98.01229 23.68574
[51,]  97.21860 23.21080
[52,]  96.52970 21.22388
[53,]  98.02872 20.26507
[54,]  97.92760 21.39398
[55,]  97.03704 20.46234
[56,]  96.92822 20.85826
[57,]  98.47076 14.08664
[58,]  98.77515 10.98880
[59,]  96.16459 21.96431
[60,]  95.97341 20.95960
[61,]  95.54920 21.47201
[62,]  95.15564 20.94652
[63,]  96.25009 22.62753
[64,]  96.06338 20.49193
[65,]  97.84218 15.79499
[66,]  97.26005 17.13186
[67,]  96.16788 19.63929
[68,]  96.13921 20.04539
[69,]  93.92534 19.67628
[70,]  98.93770 12.36708
[71,]  96.24776 16.89759
[72,]  96.03014 17.26904
[73,]  96.28511 16.66221
[74,]  96.13431 16.82926

4.12.1.2 Determining the cut-off distance

For the fixed distance weights, first, we need to determine the upper limit for distance band by using the steps below:

• Return a matrix with the indices of points belonging to the set of the k nearest neighbours of each other by using knearneigh() of spdep.

• Convert the knn object returned by knearneigh() into a neighbours list of class nb with a list of integer vectors containing neighbour region number ids by using knn2nb().

• Return the length of neighbour relationship edges by using nbdists() of spdep. The function returns in the units of the coordinates if the coordinates are projected, in km otherwise.

• Remove the list structure of the returned object by using unlist().

k1 <- knn2nb(knearneigh(coords))
k1dists <- unlist(nbdists(k1, coords, longlat = TRUE))

summary(k1dists)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  14.26   49.33   66.03   71.79   82.19  196.85 

The summary report shows that the largest first nearest neighbour distance is 196.85 km, so using this as the upper threshold gives certainty that all units will have at least one neighbour.

wm_d197 <- dnearneigh(coords, 0, 197, longlat = TRUE)
wm_d197

Neighbour list object:
Number of regions: 74 
Number of nonzero links: 828 
Percentage nonzero weights: 15.12053 
Average number of links: 11.18919 
2 disjoint connected subgraphs

Next, nb2listw() is used to convert the nb object into spatial weights object.

wm197_lw <- nb2listw(wm_d197, style = 'B')
summary(wm197_lw)

Characteristics of weights list object:
Neighbour list object:
Number of regions: 74 
Number of nonzero links: 828 
Percentage nonzero weights: 15.12053 
Average number of links: 11.18919 
2 disjoint connected subgraphs
Link number distribution:

 1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 
 3  1  2  3  3  5  3  3  7  3  4  4  6  5  5  4  1  5  4  3 
3 least connected regions:
16 57 58 with 1 link
3 most connected regions:
39 42 62 with 20 links

Weights style: B 
Weights constants summary:
   n   nn  S0   S1    S2
B 74 5476 828 1656 45160

4.12.2 Computing Adaptive distance weight matrix

One of the characteristics of fixed distance weight matrix is that more densely settled areas (usually the urban areas) tend to have more neighbours and the less densely settled areas (usually the rural counties) tend to have lesser neighbours.

Having many neighbours smoothes the neighbour relationship across more neighbours.

However, it is also possible to control the numbers of neighbours directly using k-nearest neighbours, either accepting asymmetric neighbours or imposing symmetry as shown in the code below.

knn <- knn2nb(knearneigh(coords, k=8))
knn

Neighbour list object:
Number of regions: 74 
Number of nonzero links: 592 
Percentage nonzero weights: 10.81081 
Average number of links: 8 
Non-symmetric neighbours list

Next, nb2listw() is used to convert the nb object into spatial weights object.

knn_lw <- nb2listw(knn, style = 'B')
summary(knn_lw)

Characteristics of weights list object:
Neighbour list object:
Number of regions: 74 
Number of nonzero links: 592 
Percentage nonzero weights: 10.81081 
Average number of links: 8 
Non-symmetric neighbours list
Link number distribution:

 8 
74 
74 least connected regions:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 with 8 links
74 most connected regions:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 with 8 links

Weights style: B 
Weights constants summary:
   n   nn  S0   S1    S2
B 74 5476 592 1036 19636

4.12.3 Computing GI statistics - Fixed distance

gi.fixed <- localG(Battles_2023$Incidents, wm197_lw)
gi.fixed

 [1] -2.16251076 -2.26945742 -2.36608309 -2.28060756 -1.40482297 -1.73401341
 [7] -1.78091543 -2.09517182  1.58404920  3.74739192  1.77357701  1.43327980
[13]  1.76425408  0.23868264 -0.46980525  0.66543351 -0.74168975 -1.33781214
[19] -0.51148627 -0.52707475  0.62973486  0.75874852  2.45091315 -0.79794591
[25] -0.46714759  0.78360503 -2.26619716  1.90728203 -0.19741884  0.53880930
[31] -1.05349256 -1.42802753 -0.14931467  3.57984196  2.13682831  0.42422030
[37]  1.15690841  2.06391206  2.02017424  0.24781966  2.12585948  2.66825096
[43]  0.59177302  2.14376526  1.90967745  1.33161281  1.56642737  0.66429027
[49]  3.18467135 -0.08013071  2.14013808  0.26220228 -0.29653380 -0.70655951
[55] -1.52151868 -1.55549239  1.38510282  1.35703308  1.80376535  0.80349295
[61]  1.58567210  0.37028741  1.42615140 -0.19363169  0.39364611 -0.86461442
[67] -1.28913769 -1.59126895 -1.85739540  0.40386411 -2.04421053 -1.82676914
[73] -1.82291155 -2.01994548
attr(,"internals")
               Gi      E(Gi)        V(Gi)       Z(Gi) Pr(z != E(Gi))
 [1,] 0.068229167 0.17808219 0.0025805215 -2.16251076   0.0305788284
 [2,] 0.007285974 0.09589041 0.0015242874 -2.26945742   0.0232405238
 [3,] 0.020046863 0.12328767 0.0019038942 -2.36608309   0.0179774092
 [4,] 0.006769071 0.09589041 0.0015270818 -2.28060756   0.0225716798
 [5,] 0.094095941 0.16438356 0.0025033091 -1.40482297   0.1600739255
 [6,] 0.110194304 0.20547945 0.0030195729 -1.73401341   0.0829157041
 [7,] 0.065091864 0.15068493 0.0023098862 -1.78091543   0.0749262673
 [8,] 0.091196626 0.20547945 0.0029752444 -2.09517182   0.0361557216
 [9,] 0.193083573 0.12328767 0.0019414335  1.58404920   0.1131825243
[10,] 0.289515279 0.12328767 0.0019676510  3.74739192   0.0001786828
[11,] 0.202220460 0.12328767 0.0019806822  1.77357701   0.0761331435
[12,] 0.267664828 0.19178082 0.0028030997  1.43327980   0.1517778927
[13,] 0.219305224 0.13698630 0.0021770936  1.76425408   0.0776892110
[14,] 0.091077575 0.08219178 0.0013859604  0.23868264   0.8113516776
[15,] 0.040245203 0.05479452 0.0009590683 -0.46980525   0.6384941605
[16,] 0.023995827 0.01369863 0.0002394576  0.66543351   0.5057732553
[17,] 0.090454904 0.12328767 0.0019596135 -0.74168975   0.4582753304
[18,] 0.074313409 0.13698630 0.0021946699 -1.33781214   0.1809576830
[19,] 0.139005236 0.16438356 0.0024618296 -0.51148627   0.6090105985
[20,] 0.125490196 0.15068493 0.0022849420 -0.52707475   0.5981416802
[21,] 0.058130190 0.04109589 0.0007317001  0.62973486   0.5288680718
[22,] 0.078141499 0.05479452 0.0009468162  0.75874852   0.4480030085
[23,] 0.340266667 0.20547945 0.0030244162  2.45091315   0.0142494332
[24,] 0.200730880 0.24657534 0.0033008582 -0.79794591   0.4249018788
[25,] 0.167275574 0.19178082 0.0027517563 -0.46714759   0.6403942853
[26,] 0.303858068 0.26027397 0.0030935821  0.78360503   0.4332719050
[27,] 0.062418386 0.17808219 0.0026049511 -2.26619716   0.0234393145
[28,] 0.355729167 0.24657534 0.0032752773  1.90728203   0.0564840763
[29,] 0.061812467 0.06849315 0.0011451559 -0.19741884   0.8434997910
[30,] 0.146966527 0.12328767 0.0019313068  0.53880930   0.5900184491
[31,] 0.043455497 0.08219178 0.0013519884 -1.05349256   0.2921153017
[32,] 0.030184751 0.08219178 0.0013263280 -1.42802753   0.1532839350
[33,] 0.076701571 0.08219178 0.0013519884 -0.14931467   0.8813053367
[34,] 0.363684489 0.17808219 0.0026880602  3.57984196   0.0003438021
[35,] 0.322002635 0.20547945 0.0029736197  2.13682831   0.0326119589
[36,] 0.171367177 0.15068493 0.0023769086  0.42422030   0.6714051589
[37,] 0.252951981 0.19178082 0.0027957315  1.15690841   0.2473097820
[38,] 0.232058669 0.13698630 0.0021219065  2.06391206   0.0390260550
[39,] 0.396593674 0.27397260 0.0036842794  2.02017424   0.0433653170
[40,] 0.047619048 0.04109589 0.0006928579  0.24781966   0.8042739429
[41,] 0.387329591 0.26027397 0.0035720591  2.12585948   0.0335149615
[42,] 0.435444414 0.27397260 0.0036621834  2.66825096   0.0076247282
[43,] 0.134509081 0.10958904 0.0017733202  0.59177302   0.5540025915
[44,] 0.370217451 0.24657534 0.0033264297  2.14376526   0.0320517005
[45,] 0.132483082 0.06849315 0.0011228022  1.90967745   0.0561747560
[46,] 0.113924051 0.06849315 0.0011639833  1.33161281   0.1829874538
[47,] 0.307425214 0.21917808 0.0031738082  1.56642737   0.1172486006
[48,] 0.137802607 0.10958904 0.0018038490  0.66429027   0.5065045498
[49,] 0.339932274 0.17808219 0.0025828347  3.18467135   0.0014491849
[50,] 0.092809365 0.09589041 0.0014784223 -0.08013071   0.9361332989
[51,] 0.287356322 0.17808219 0.0026070606  2.14013808   0.0323436091
[52,] 0.261594581 0.24657534 0.0032811258  0.26220228   0.7931654986
[53,] 0.071372753 0.08219178 0.0013311533 -0.29653380   0.7668224621
[54,] 0.080156658 0.10958904 0.0017352153 -0.70655951   0.4798402603
[55,] 0.123498695 0.20547945 0.0029031487 -1.52151868   0.1281297272
[56,] 0.131920530 0.21917808 0.0031468082 -1.55549239   0.1198288467
[57,] 0.035637728 0.01369863 0.0002508843  1.38510282   0.1660210309
[58,] 0.034807642 0.01369863 0.0002419663  1.35703308   0.1747706992
[59,] 0.366230366 0.26027397 0.0034505972  1.80376535   0.0712681006
[60,] 0.292600313 0.24657534 0.0032811258  0.80349295   0.4216898674
[61,] 0.354370214 0.26027397 0.0035214197  1.58567210   0.1128137120
[62,] 0.295910393 0.27397260 0.0035100061  0.37028741   0.7111683500
[63,] 0.299541655 0.21917808 0.0031753179  1.42615140   0.1538246447
[64,] 0.222019781 0.23287671 0.0031438461 -0.19363169   0.8464642818
[65,] 0.066967845 0.05479452 0.0009563271  0.39364611   0.6938423348
[66,] 0.108660999 0.15068493 0.0023623728 -0.86461442   0.3872504579
[67,] 0.124184712 0.19178082 0.0027494435 -1.28913769   0.1973502243
[68,] 0.131770833 0.21917808 0.0030172252 -1.59126895   0.1115490605
[69,] 0.041992697 0.12328767 0.0019156610 -1.85739540   0.0632549198
[70,] 0.036378335 0.02739726 0.0004945225  0.40386411   0.6863126506
[71,] 0.063607925 0.16438356 0.0024302999 -2.04421053   0.0409327537
[72,] 0.096329081 0.19178082 0.0027302372 -1.82676914   0.0677344872
[73,] 0.085438916 0.17808219 0.0025828347 -1.82291155   0.0683167878
[74,] 0.065070276 0.16438356 0.0024173270 -2.01994548   0.0433890436
attr(,"cluster")
 [1] Low  Low  Low  Low  Low  High Low  Low  Low  Low  High Low  High High High
[16] Low  Low  High Low  Low  High High High Low  Low  High Low  Low  Low  Low 
[31] Low  Low  Low  High Low  High Low  Low  High Low  High High Low  High Low 
[46] High High High Low  High Low  Low  Low  Low  Low  High High Low  Low  Low 
[61] High Low  High Low  High High Low  Low  Low  High Low  Low  Low  Low 
Levels: Low High
attr(,"gstari")
[1] FALSE
attr(,"call")
localG(x = Battles_2023$Incidents, listw = wm197_lw)
attr(,"class")
[1] "localG"

Next, we will join the Gi values to their corresponding sf data frame by using the code below.

Battles_2023.gi <- cbind(Battles_2023, as.matrix(gi.fixed)) %>%
  rename(gstat_fixed = as.matrix.gi.fixed.)

Battles_2023.gi

Simple feature collection with 74 features and 11 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 99.66532 ymax: 28.54554
Geodetic CRS:  WGS 84
First 10 features:
           DT   DT_PCODE    DT_MMR PCode_V year event_type Incidents Fatalities
1    Hinthada MMR017D002  ဟင်္သာတခရိုင်     9.4 2023    Battles         4          3
2     Labutta MMR017D004  လပွတ္တာခရိုင်     9.4 2023    Battles         1          1
3     Pathein MMR017D001    ပုသိမ်ခရိုင်     9.4 2023    Battles         3          1
4      Pyapon MMR017D006   ဖျာပုံခရိုင်     9.4 2023    Battles         3          2
5        Bago MMR007D001    ပဲခူးခရိုင်     9.4 2023    Battles        50        270
6     Taungoo MMR007D002  တောင်ငူခရိုင်     9.4 2023    Battles        87        493
7        Pyay MMR008D001    ပြည်ခရိုင်     9.4 2023    Battles        34         59
8  Thayarwady MMR008D002 သာယာဝတီခရိုင်     9.4 2023    Battles        50         93
9       Falam MMR004D001   ဖလမ်းခရိုင်     9.4 2023    Battles        27         89
10      Hakha MMR004D003 ဟားခါးခရိုင်     9.4 2023    Battles        48        210
   Z.Incidents lag_Incidents gstat_fixed                       geometry
1  -0.80534765      18.20000   -2.162511 MULTIPOLYGON (((95.12637 18...
2  -0.85573862       3.00000   -2.269457 MULTIPOLYGON (((95.04462 15...
3  -0.82214464       2.50000   -2.366083 MULTIPOLYGON (((94.27572 15...
4  -0.82214464       3.00000   -2.280608 MULTIPOLYGON (((95.20798 15...
5  -0.03268604      40.66667   -1.404823 MULTIPOLYGON (((95.90674 18...
6   0.58880265      29.00000   -1.734013 MULTIPOLYGON (((96.17964 19...
7  -0.30143790      31.40000   -1.780915 MULTIPOLYGON (((95.70458 19...
8  -0.03268604      35.60000   -2.095172 MULTIPOLYGON (((95.85173 18...
9  -0.41901684      53.00000    1.584049 MULTIPOLYGON (((93.36931 24...
10 -0.06628002      62.00000    3.747392 MULTIPOLYGON (((93.35213 23...

Note

The codes above performs three tasks.

• First, it converts the output vector (i.e. gi.fixed) into r matrix object by using as.matrix().

• Next, cbind() is used to join Battles_2023 and gi.fixed matrix to produce a new SpatialPolygonDataFrame called Battles_2023.gi.

• Lastly, the field name of the gi values is renamed to gstat_fixed by using rename().

4.12.3.1 Mapping Gi values with Fixed distance weights

The code below plots the Gi values derived using fixed distance weight matrix, for event type==Battles in 2023.

Gimap <-tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +
  
  tm_shape(Battles_2023.gi) +
  tm_fill(col = "gstat_fixed", 
          style = "pretty",
          palette="-RdBu",
          title = "Fixed Distance\nlocal Gi") +
  tm_borders(alpha = 0.5)

Gimap

4.12.4 Computing GI statistics - Adaptive distance

The code below is used to compute the Gi values for Incidents of Battles in 2023 by using an adaptive distance weight matrix (i.e knb_lw).

gi.adaptive <- localG(Battles_2023$Incidents, knn_lw)
Battles_2023.gi <- cbind(Battles_2023, as.matrix(gi.adaptive)) %>%
  rename(gstat_adaptive = as.matrix.gi.adaptive.)

datatable(Battles_2023.gi)

4.12.4.1 Mapping Gi values with Adaptive distance weights

The code below plots the Gi values derived using adaptive distance weight matrix for event type == Battles in 2023.

Gimap <- tm_shape(mmr_shp_mimu_2) +
  tm_borders() +  # Draws borders for all regions
  tm_fill(col = "white", alpha = 0.5, title = "Background") +
  
  tm_shape(Battles_2023.gi) + 
  tm_fill(col = "gstat_adaptive", 
          style = "pretty", 
          palette="-RdBu", 
          title = "Adaptive Distance\nlocal Gi") + 
  tm_borders(alpha = 0.5)

Gimap

Parameters and Output to be exposed to Shiny

In terms of enabling user inputs, I will expose it to accept user inputs:

• Year

• Event type

• Data Classification type, and

• Number of clusters for the Adaptive weight matrix

4.13 UI design - Part 4

Both the choropleth map and data table above will be used to populate our app’s Hot & Cold Spot Analysis module.

The motivation here is to provide users with a visualisation of statistically significant areas of event types in order to help users understand how hot and cold zones may have developed and/or changed over time.

The below is a visual of the prototype.

Functionality and interactivity

• Users will be able to select the specific year, event type, the data classification method and the number of clusters for the adaptive weight matrix to populate the map and data table.

Additional considerations

Here, I have chosen to implement a single point for users to filter and select, enabling users to see the map and refer to data table with the corresponding Gi statistics.

Additionally, I will first implement adaptive distance weights, instead of having both fixed and adaptive. This is due the the added complexity of the shiny code required to accomodate the selection of either adaptive or distance weights within the same reactive statement.

For example, as shown above in the calculation for fixed distance weights, we will first need to determine the upper limit of the distance band in order to set such that all units will have at least one neighbour.

When considering that user inputs selections will require a recalculation of the dataset to calculate the GI values reactively, this may add complexity to the code and the computation required to plot the map.

4.14 References

Main reference: Kam, T.S. (2024). Local Measures of Spatial Autocorrelation

Imran’s Data Science website Imran’s Data Science website Exercises in R R Exercise 1 - Analysing PISA’s education survey data /R-ex/R-Ex1/R_Ex1.html R Exercise 2 - 5 Exploratory Data Analyses on Pisa Data /R-ex/R-Ex2/R_Ex2.html R Exercise 3 - DataVis Makeover of another student’s work /R-ex/R-Ex3/R_Ex3.html R Exercise 4 - Be Weatherwise or Otherwise /R-ex/R-Ex4/R_Ex4.html R Exercise 5 - Visual Analytics exercise on Armed conflicts - Initial Draft /R-ex/R-Ex5/R_Ex5.html Exercises in Python Python Exercise 1 /Python-ex/Python-Ex1/Python_Ex1.html Python Exercise 2 Python Exercise 3 Write-ups Write-up 1 /Write-ex/Write-Ex1/Write_Ex1.html Write-up 2 /Write-ex/Write-Ex2/Write_Ex2.html Write-up 3 /Write-ex/Write-Ex3/Write_Ex3.html /about.html https://github.com/imranmi/imran-s-data-sc https://www.linkedin.com/in/imran-ibrahim-116a236b

Imran’s Data Science website - Visual Analytics exercise on Armed conflicts - Initial Draft Imran’s Data Science website - Visual Analytics exercise on Armed conflicts - Initial Draft Imran’s Data Science website - Visual Analytics exercise on Armed conflicts - Initial Draft Imran’s Data Science website